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Abstract 

The chaperonin GroEL-GroES, a machine which helps some proteins to fold, cycles through a number 
of allosteric states, the T state, with high affinity for substrate proteins (SPs), the ATP-bound R state, 
and the R" (GroEL — ADP — GroES) complex. Structures are known for each of these states. Here, 
we use a self-organized polymer (SOP) model for the GroEL allosteric states and a general structure- 
based technique to simulate the dynamics of allosteric transitions in two subunits of GroEL and the 
heptamer. The T ^ R transition, in which the apical domains undergo counter-clockwise motion, 
is mediated by a multiple salt-bridge switch mechanism, in which a series of salt-bridges break and 
form. The initial event in the R — > R" transition, during which GroEL rotates clockwise, involves a 
spectacular outside-in movement of helices K and L that results in K80-D359 salt-bridge formation. 
In both the transitions there is considerable heterogeneity in the transition pathways. The transition 
state ensembles (TSEs) connecting the T, R, and R" states are broad with the the TSE for the T ^ R 
transition being more plastic than the R R" TSE. The results suggest that GroEL functions as a 
force-transmitting device in which forces of about (5-30) pN may act on the SP during the reaction 
cycle. 
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INTRODUCTION 



The hallmark of allostery in biomolecules is the conformational changes at distances far 
from the sites at which ligands bind AUosteric transitions arise from long-range spa- 

tial correlations between a network of residues. Signaling between the residues likely results 
in conformational transitions. The potential link between large scale allosteric transitions and 

n 

function is most vividly illustrated in biological nanomachines [3|. For example, during tran- 
scription DNA polymerases undergo a transition from open to closed state that is triggered by 
dNTP binding [4]. Although allosteric transitions were originally associated with multi-subunit 
assemblies, it is now believed that a network of residues encode the dynamics of monomeric 
globular proteins [sj]. 

Computational methods have been used to determine the network of spatially correlated 
residues (wiring diagram) that trigger the functionally relevant conformational changes. Re- 
cently, static methods, either sequence ^,[7,0] or structure-based Q, [llj] have been proposed 
to predict the allosteric wiring diagram. However, in order to fully understand the role of al- 



lostery it is important to dynamical 



from one state to another 



0, y, E, 



monitor the structural changes that occur in the transition 



10|. Here, we propose a method for determining the al- 



losteric mechanism in biological systems with applications to dynamics of such processes in the 
chaperonin GroEL, an ATP-fueled nanomachine, which facilitates folding of proteins (SPs) that 

nn 

are otherwise destined to aggreg ate 116|,117[. 

GroEL consists of two heptameric rings, stacked back-to-back. Substrate proteins are cap- 
tured by GroEL in the T state (Fig. 1) while ATP-binding triggers a transition to the R state. 
The T and R structures show that the equilibrium T R transition results in large (nearly 
rigid body) movements of the apical (A) domain (Fig. 1), and somewhat smaller changes in 
the conformations of the intermediate (I) domain. Binding of the co-chaperonin GroES requires 
dramatic movements in the A domains which double the volume of the central cavity. Compar- 
ison of the structures of the T, R, and the R" {GroEL — {ADP)j — GroES) indicates that the 
equatorial (E) domain, which serves as an anchor jl6l|, undergoes comparatively fewer structural 
changes. Although structural and mutational studies 18|, ll9|, l20| have identified many important 
residues that affect GroEL function, only few studies have explored the dynamics of allosteric 
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transitions between the various states |2l|, |22|, |23 |. 



Here, we use the self-organized polymer (SOP) model of GroEL and a novel technique 
(Methods) to monitor the order of events in the T ^ R, R ^ R", and T R" transitions. 
By simulating the dynamics of ligand-induced conformational changes in the heptamer and two 
subunits we have obtained an unprecedented view of the key interactions that drive the various 
allosteric transitions. The dynamics of transition between the states are achieved (Methods) 
under the assumption that the rate of conformational changes is slower than the rate at which 
ligand-binding induced strain propagates. Because of the simplicity of the SOP model we have 
been able to generate multiple trajectories to resolve the key events in the allosteric transitions. 
We make a number of predictions including the identification of a multiple salt-bridge switch 
mechanism in the T ^ R transition, and the occurrence of dramatic movement of helices K 
and L in the R R" transition. The structures of the transition state ensembles that connect 
the various end points show considerable variability mostly localized in the A domain. 

RESULTS and DISCUSSION 

Heptamer dynamics show that the A domains rotate counter-clockwise in the T 
R transition and clockwise in the R R" transition: In order to probe the global motions 
in the various stages of GroEL allostery we simulated the entire heptamer (Methods). The 
dynamics of the T R transition, monitored using the time-dependent changes in the angles 
a, /3, and 7 (see the caption in Fig. 2 for definitions), that measure the relative orientations 
of the subunits, show (Fig. 2-A) that the A domains twist in a counter-clockwise manner in 



agreement with experiment [2^. The net changes in the angles in the R — R" transition, which 
occurs in a clockwise direction (Fig. 2-B), is greater than in the T R transition. As a result 
the global T R" transition results in a net ~ 110° rotation of the A domains. Surprisingly, 
there are large variations in the range of angles explored by the individual subunits during the 
T ^ R ^ R" transitions. There are many more inter-subunit contacts in the E domain than 
in the A domain, thus permitting each A domain to move more independently of one another. 
Fig. 2 shows that the dynamics of each subunit is distinct despite the inference, from the end 
states alone, that the overall motion occurs without significant change in the root mean square 



deviation (RMSD) of the individual domains. The time- dependent changes in the angles a, 
/?, and 7 from one subunit to another are indicative of an inherent dynamic asymmetry in the 
individual subunits that has been noted in static structures 
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26(1 . As in the T ^ R transition, 



there is considerable dispersion in the time-dependent changes in a, (3, and 7 of the individual 
subunits (Fig. 2-B) during the R R" transition. 

The clockwise rotation of apical domain alters the nature of lining of the SP binding sites 
(domain color-coded in magenta in Fig. 1). The dynamic changes in the 7 angle (Fig. 2) 
associated with the hinge motion of the I domain that is perpendicular to the A domain lead 
to an expansion of the overall volume of the heptamer ring. More significant conformational 
changes, that lead to doubling of the volume of the cavity, take place in the R R" transition. 
The apical domain is erected, so that the SP binding sites are oriented upwards providing 
binding interfaces for GroES. Some residues, notably 357-361, which are completely exposed 
on the exterior surface in the T state move to the interior surface during the T R ^ R" 
transitions. 



Global T R and R — > R" transitions follow two-state kinetics: Time-dependent 
changes in RMSD with respect to a reference state (T, R, or R"), from which a specific 
allosteric transition commences (Fig. 3), differ from molecule-to-molecule, which reflects 
the heterogeneity of the underlying dynamics (Fig. 3-A). Examination of the RMSD for 
a particular trajectory in the transition region (Fig. 3-A inset) shows that the molecule 
undergoes multiple passages through the transition state (TS). Assuming that RMSD is a 
reasonable reaction coordinate, we find that GroEL spends a substantial fraction of time 
(measured with respect to the first passage time) in the TS region during the T — > i? 
transition. By averaging over 50 individual trajectories we find that the ensemble average of 
the time-dependence of RMSD for both the T ^ R and R R" transitions follow single 
exponential kinetics. From this observation we conclude that, despite a broad transition 
region, the allosteric transitions can be approximately described by a two-state model. From 
the exponential fits of the global RMSD changes we find that the relaxation times for the 
two transitions are tt^r ~ 25fis and tr^rh ~ 71/is. The larger time constant for the 
R — > R" transition compared to the T ^ R transition is due to the substantially greater 

4 



rearrangement of the structure in the former. Unhke the global dynamics characterizing the 
overall motion of GroEL, the local dynamics describing the formation and rupture of key 
interactions that drive GroEL allostery cannot be described using two-state kinetics (see below). 



T — > i? transition is triggered by downward tilt of helices F and M in the I-domain 
followed by a multiple salt-bridge switching mechanism: Several residues in helices F 
(141-151) and M (386-409) in the I domain interact with the nucleotide-binding sites in the E 
domain thus creating a tight nucleotide binding pocket. The favorable interactions are enabled 

ainding 
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by the F, M helices tilting by about 15° that results in the closing of the nucleotide 
sites. A number of residues around the nucleotide binding pocket are highly conserved 
Since the T —>■ R transition involves the formation and breakage of intra- and inter- subunit 
contacts we simulated two interacting subunits so as to dissect the order of events. 

(i) The ATP-binding-induced downward tilt of the F, M helices is the earliest event [22] that 
accompanies the subsequent spectacular movement of GroEL. The changes in the angles F and 
M helices make with respect to their orientations in the T state occur in concert (Fig. 3-C). At 
the end of the R R" transition the helices have tilted on average by about 25° in all (Fig. 
3-C) . There are variations in the extent of the tilt depending on the molecule (see inset in Fig. 
3-G). Upon the downward tilt of the F and M helices, the entrance to the ATP binding pocket 
narrows as evidenced by the rapid decrease in the distance between P33 and N153 (Fig. 4). 
A conserved residue, P33 contacts ATP in the T state and ADP in the R" structure and is 



involved in allostery 



231]. The contact number of N153 increases substantially as a result of loss 



in accessible surface area during the R R" transition 5 



. In the T state E386, located at the 



tip of M helix, forms inter-subunit salt-bridges with R284, R285, and R197. In the transition to 
the R state these salt-bridges are disrupted and the formation of a new intra-subunit salt-bridge 
with K80 takes place simultaneously (see the middle panel in Fig. 4). The dynamics show 
that the tilting of M helix must precede the formation of inter-subunit salt-bridge between the 
charged residues E386 with K80. 

(ii) The rupture of the intra-subunit salt-bridge D83-K327 occurs nearly simultaneously with 
the disruption of the E386-R197 inter-subunit interaction. The distance between the Ca atoms 
of D83 and K327 of around 8.5 A in the T state slowly increases to the equilibrium distance 
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of ~ 13 A in the R state with relaxation time r ^ 100 j2s (the top panel in Fig. 4). The 
establishment of K80-E386 salt-bridge occurs around the same time as the rupture of R197- 
E386 interaction. In the T — i? formation a network of salt-bridges are broken and new ones 
formed (see below). At the residue level, the reversible formation and breaking of D83-K327 
salt-b ridge, i n concert with the inter-subunit salt-bridge switch associated with E386 |24] and 



E257 
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30l | , are among the most significant events that dominate the T R transition. 



Remarkably, the coordinated global motion is orchestrated by a multiple salt-bridge switching 
mechanism. The movement of the A domain results in the dispersion of the SP binding sites 
(Fig. 1) and also leads to the rupture of the E257-R268 inter-subunit salt-bridge. The kinetics 
of breakage of the E257-R268 salt-bridge are distinctly non-exponential (the last panel in Fig. 
4). It is very likely that the dislocated SP binding sites maintain their stability through the 
inter-subunit salt-bridge formation between the apical domain residues. To maintain the stable 
configuration in the R state, E257 engages in salt-bridge formation with positively charged 
residues that are initially buried at the interface of inter-apical domain in the T state. Three 
positively charged residue at the interface of the apical domain in the R state, namely, K245, 
K321, and R322 are the potential candidates for the salt-bridge with E257. During the T ^ R 
transitions E257 interact partially with K245, K321, and R322 as evidenced by the decrease in 
their distances (the last panel in the middle column of Fig. 4). The distance between E409-R501 
salt-bridge, involving residues that connect and hold the I and E domains, remains intact at a 
distance ~ 10 A throughout the whole allosteric transitions. This salt-bridge and two others 
(E408-K498 and E409-K498) might be important for enhancing positive intra-ring cooperativity 
and for stability of the chaperonins. Indeed, mutations at sites E409 and R501 alter the stability 
of the various allosteric states jsi]. In summary, we find a coordinated dynamic changes in the 
network of salt-bridges are linked in the T ^ R transition. 

The order of events, described above, are not followed in all the trajectories. Each 
molecule follows somewhat different pathway during the allosteric transitions which is indi- 
cated by the considerable dispersion in the dynamics. However, when the time traces are 
averaged over a large enough sample, the global kinetics can be analyzed using a two state model. 



R — > R" transition involves a spectacular outside-in movement of K and L helices 
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accompanied by inter-domain salt-bridge formation K80-D359: The dynamics of the 
irreversible R — > R" transition is propelled by substantial movements in the A domain helices 

K and L that drive the dramatic conformational change in GroEL resulting in doubling of the 
volume of the cavity The dynamics of the R R" transition also occur in stages. 

(i) Upon ATP hydrolysis the F, M helices rapidly tilt by an additional 10° (Fig. 3-C). Nearly 
simultaneously there is a small reduction in P33-N153 distance (7 A— > 5 A) (see top panel in 
Fig. 5). These relatively small changes are the initial events in the R R" transition. 

(ii) In the next step, the A domain undergoes significant conformational changes that arc 
most vividly captured by the outside- in concerted movement of helices K and L. Helices K and 
L, that tilt by about 30" during the T — > transition, further rotate by an additional 40" when 
the R — > R" transition occurs (Fig. 3-D). In the process, a number of largely polar and charged 
residues that are exposed to the exterior in the T state line the inside of the cavity in the R" 
state. The outside-in motion of K and L helices leads to an inter-domain salt-bridge K80-D359 
whose Ca distance changes rapidly from about 40 A in the R state to about 14 A in the R" 
(Fig. 5). 

The wing of the apical domain that protrudes outside the GroEL ring in the R state moves 
inside the cylinder. The outside-in motion facilitates the K80-D359 salt-bridge formation which 
in turn orients the position of the wing. The orientation of the apical domain's wing inside the 
cyhnder exerts a substantial strain (data not shown) on the GroEL structure. To relieve the 
strain, the apical domain is forced to undergo a dramatic 90" clockwise rotation and 40" upward 
movement with respect to the R state. As a result, the SP binding sites (H, 1 hehces) are 
oriented in the upward direction. Before the strain-induced alterations are possible the distance 
between K80 and D359 decreases drastically from that in R state (middle panel in Fig. 5) . The 
clockwise motion of the apical domain occurs only after the formation of salt-bridge between K80 
and D359. On the time scale during which K80-D359 salt-bridge forms, the rupture kinetics of 
several inter-apical domain salt-bridges involving residues K245, E257, R268, K321, and R322, 
follow complex kinetics (Fig. 5). Formation of contact between 1305 and A260 (a binding site 
for substrate proteins), and inter-subunit residue pair located at the interface of two adjacent 
apical domains in the R" state, occurs extremely slowly compared to others. The non-monotonic 
and lag-phase kinetics observed in the rupture and formation of a number of contacts suggests 
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that intermediate states must exist in the pathways connecting the R and R" states . 

The clockwise rotation of apical domain, that is triggered by a network of salt-bridges as 
well as interactions between hydrophobic residues at the interface of subunits, orients it in 
the upward direction so as to permit the binding of the mobile loop of GroES. Hydrophobic 
interactions between SP binding sites and GroES drive the R —>■ R" transition. The hydrophilic 
residues, that are hidden on the side of apical domain in the T or the R state, are now exposed 
to form an interior surface of the GroEL (see the residue colored in yellow on the A domain 
in Fig. 1). The E409-R501 salt-bridge formed between I and A domains close to the 7 — Pj 

binding site is maintained throughout the allosteric transitions including in the transition state 

31|. 

Transition state ensembles (TSEs) are broad: The structures of the TSEs connecting 
the T, i?, and R" states are obtained using RMSD as a surrogate reaction coordinate. We assume 
that, for a particular trajectory, the TS location is reached when 5^ = \{RMSD /T){tTs) -~ 
[RMSD / R){tTs)\ < "^c where = 0.2 A, and Its is the time at which 5^ < r^. Letting the value 
of RMSD at the TS be = l/2x\{RMSD/T){tTs) + {RMSD/R){tTs)\ the distributions P(A*) 
ioT T —>■ R and R R" transitions are broad (see Fig. S3 in the Supplementary Information). 
If is normalized by the RMSD between the two end point structures to produce a Tanford 
/9-like parameter (see caption to Fig. 6 for definition), we find that the width of the TSE for 
the R R" is less than the T R transition (Fig. 6-A). The mean values of for the two 
transitions show that the most probable TS is located close to the R states in both T ^ R and 
R R" transitions. 

Disorder in the TSE structures (Fig. 6) is largely localized in the A domain which shows that 
the substructures in this domain partially unfold as the barrier crossings occur. By comparison 
the E domain remains more or less structurally intact even at the transition state which 
suggests that the relative immobility of this domain is crucial to the function of this biological 
namomachine [3]. The dispersions in the TSE are also reflected in the heterogeneity of the dis- 
tances between various salt-bridges in the transition states. The values of the contact distances, 
in the T R transition among the residues involved in the salt-bridge switching between K80, 
R197, and E386 at the TS has a very broad distribution (Fig. 6-B) which also shows that 
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the R197-E386 is at least partially disrupted in the TS and the K80-E386 is partially formed 



32] 



CONCLUSIONS 



Based on the observation that the conformational changes in molecular nanomachines are 
much slower than the rate of ligand binding-induced strain propagation, we have developed a 
structure-based method to probe the allosteric transitions in biological molecules. In our method, 
transitions between multiple states are probed using the state-dependent energy functions and 
Brownian dynamics simulations. Applications to allosteric transitions in GroEL have produced 
a number of new predictions that can be experimentally tested. The global dynamics in the 
T R transition reveals that the A domain rotates in a counter-clockwise manner whereas they 
rotate in a clockwise direction during the R R" transition. Although this observation was 
previously made based on the end point structures (T, i?, and R") alone, the transition kinetics 
show substantial dynamic heterogeneity in the dynamics of individual subunits. A key finding 
of the present work is that the transitions occur by a multiple coordinated switch between a 
network of salt-bridges. Our findings suggest a series of mutational studies that can examine 
the link between disruption of the salt-bridge interactions and the GroEL function. The most 
dramatic outside-in movement, the rearrangement of helices K and L of the A domain, occurs 
largely in the R — * R" transition and results in the inter-subunit K80-D359 salt-bridge formation. 
The TSEs are broad with the T R TSE being more plastic than the one connecting R and 
R" states. In both the transitions most of the conformational changes occur in the A domain 
with the E domain serving as a largely structurally unaltered static base that is needed for force 
transmission jl^. 

The unprecedented picture of the dynamics of allostery in GroEL presented here suggests 
that, like other ATP-consuming nanomachines, GroEL may be viewed as a force-generating 
device. The combination of large-scale (nearly ~ 50 A) dispersion in the SP binding sites 
and the multivalent binding of SP must axiomatically lead to considerable strain on the SP. 
Mechanochemical considerations enable us to make a rough estimate of the equilibrium force, /, 
acting on the SP resulting from the observed allosteric movements from / ~ — ^^'^ — - where 
the ATP hydrolysis free energy AG^^^{= AG° + i?T log ^^^p^pp) at physiological condition 



{[ATP] = SOOuM, [ADP] = SOfxM, [Pi] = 200/iM, T = 320K, and AG" = -7.3kcal/mol) 
is ~ —12kcal / mol , and Ad ~ (10 — 50) A. Assuming an efficiency of about 0.5, we estimate 
/ ~ 5 — 30 pN which is large enough to partially or fully unfold many proteins. 



METHODS 



Using the aval 



able structures in the 



Energy function 
lOEL (T state) [33^207 

SOP Hamiltonian 3, [l3] of the states (X = T, R, Ft!') of GroEL as 
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^ 47rer,-* 



i*j*G{salt bridge} 

The ffist term, which accounts for chain connectivity is represented using finite extensible non- 
linear elastic (FENE) potential js^], with parameters /c = 20 kcal/ [mol-k.^), Rq = 2 A where 
Tj^j+i is the distance between neighboring interaction centers i and i + r° -_,_]^(X) is the distance 
in state X. The Lennard- Jones potential (second term) accounts for interactions that stabilize a 
particular allosteric state. Native contact exists if the distance between i and j less than Rc = 8 
A in state X for > 2. Hi and j sites are in contact in the native state, A^ = 1, otherwise 

Aij = 0. To ensure the non-crossing of the chain, we used a G*'^ power potential in the third and 
the fifth terms and set cr = 3.8 A, which is the — distance. We used eh = 2 kcal/mol if 
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the residues are in contact and e; = 1 kcal/mol for non-native pairs. 

The fourth and the fifth terms in Eq{T] are for interaction of residues with ATP {R state) or 
ADP {R" state). The atomic coordinates of ATP (ADP) are taken from the R {R") structure 
without coarse-graining. The functional form for residue- ATP (or ADP) interaction is the same 
as for residue-residue interactions with e^'^^ = 0.2kcal/mol, ef^^ = O.lkcal/mol. We used a 
small value of e^^^ (or ef^^) because the coordinates of all the heavy atoms of ATP and ADP 
are explicitly used as interaction sites. The distance between the i^^ residue and the j^^ atom in 
ATP (or ADP) is aij. We used the screened electrostatic potential where = 2.4 A, e = lOeo, 
and qiq2 = — to account for the favorable salt-bridge interactions which are state-independent 
(EqUD. 

Inducing allosteric transitions: The T —>■ R allosteric transition of GroEL is simulated by 
integrating the equations of motion with the force arising from H{{fi}\R). However, the ensem- 
ble of initial structures were generated using the T-state H{{fi}\T). The Brownian dynamics 



algorithm 
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391] determines the configuration of GroEL at time t as follows, 

{{} f,{t + h)=fi{t) + !^{Fi{t\T) + fi{t)) (o<t<r) 

(u) f,{t + h*) = r,{t) + ^{F,{t\T ^R) + f,(t)) (r < t < r + Nth*) (2) 

{ill) nit + h) = nit) + \{Fi{t\R) + f,{t)) {t>t* + Nth*) 

where Fi{t\X) = —V ff^H {{n}\X) (X = T, i? or T — > R), the Newtonian force acting on a 
residue i, and Ti{t) is a random force on i^^ residue that has a white noise spectrum that satisfies 

{fi{t)-fj{t + nh)) = ^6o,Ji,j, where 6o n is the Kronecker delta function and n = 0,l,2, 

As long as the fluctuation-dissipation theorem is satisfied it can be shown that our procedure 
for switching the Hamiltonian from T to R will lead to the correct Boltzmann distribution for 
the R state at long times [l4| . 

The algorithm in Eql2] is implemented in three steps: (i) During the time interval < t < t* 
an ensemble of T-state conformations is generated; (ii) The energy function is switched from 
H{{fi}\T) to H{{fi}\R) symbolized by H{{ri}\T R) in the duration t* < t < t* + Nth*. 
If Nt = our method is similar to one in 9\; (iii) A dynamic trajectory under H{{n}\R) is 
generated for t > t*. The assumption in our method is that the rate of conformational change 
in biomolecules is smaller than the rate at which a locally applied strain (due to ligand binding) 
propagates. As a result, the Hamiltonian switch should not be instantaneous [Nt ^ 0). Using 
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a non-zero value of A''^ (second step in Eql2]) not only ensures that there is a lag time between 
ligand binding and the associated response but also eliminates computational instabilities in the 
distances between certain residues that change dramatically during the transition. The "loading" 
rate can be altered by varying Nt and hence even non-equilibrium ligand-induced transitions can 
be simulated. Additional details of the simulations are given in the Supplementary Information. 

Additional details: When the transition, say T to R occurs, the native contact distances 
{r°-{T) and r°-{R)) in the two states are different. If r°j{T) -C r°j{R) then, upon the Hamiltonian 
switch, a large repulsive force arises from the the pair i and j that is initially in the equilibrated 
T-state conformation at a distance, r^j r°j{T)). To eliminate the purely computational 
problem we made a gradual transition of the parameters using a linear interpolation procedure. 
For the T ^ R transition we used, r° (T R) = _ rjj_^g al^iT i?) values 

are similarly defined with a°j{T) = but a°j{R or R") ^ 0. Accordingly, the corresponding 
Hamiltonian H{{fi}\T — > R) for stage (ii) is defined by inserting r° (T R), a°j(T R) in 
EqO We set K = 100 and increased k from to K every 6 = 100 integration time step, which 
leads to Nt = K X 6. In addition, we choose the integration time step h* = (0.001 — 0.01)tl 
during the Hamiltonian switch. Thus, the Hamiltonian switch is smoothly made within (5 — 50) 
ns without causing any computational instability. The technical problem arises only for certain 
pairs of residues for which r°-{T) <^ r°-{R). 
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FIGURE CAPTIONS 



Figure 1 : GroEL structure. The columns from left to right show T, i?, and R" structures of 
GroEL structures. The top view is given in first row (see Fig. SI in Supplementary Information 
for a side view) and the second row displays the side view of a single subunit. The white ball 
represents D359. The hehces that most directly influence the allosteric transitions transitions 
are labeled. 

Figure 2 : GroEL dynamics monitored using various angles. A. T — i? transition dynamics 
for the heptamer monitored using angles, a, 7. An angle 9 {— a, (3) is defined by cos^(i) = 

— * 

ue{^)-U0{t)/\ue{^) I \ue{t) \ . For a, we obtain Ua{t) by projecting the vector {r2m{i) (t) = -^236(1) (t) - 
Rcm) between the center of mass (Rcm) and residue 236 on i*'* subunit (-R236(i)(^)) onto the plane 
perpendicular to the principal axis (ep) of the heptamer, i.e.,. Uait) = 'r236(i)(^) — (^236(i)(^)"ep)ep. 
The angle between H helices (residue 231—242) of i*^ subunit at times t = and t using 
the vector, i?23i(i)(^) — -R242(i)(^) is (3- The sign of the angles {a and (3) is determined using 
sgn[{u{0) X u{t)) • ep], which is (+) for counter-clockwise and (-) for clockwise rotation. 7 
measures the perpendicular motion of apical domain with respect to the hinge (residue 377). 
We defined u^{t) = i?236(i)(^) ~ Rzn{i)if) at each subunit i, and 7(t) = 90" — cos^"*^ {u^ • ep). On 
the right three panels we plot the time dependence of ct, and 7 for each subunit in different 
color. The black line represents the average of 21 (=3x7) values of each angle calculated from 
three trajectories of 7 subunits. B. Same as in A except for the R — > R" transition. 

Figure 3 : RMSD as a function of time. A. Time-dependence of RMSD of a few individual 
molecules are shown for T ^ R transition. Solid (dashed) lines are for RMSD/T (RMSD/i?) 
(RMSD calculated with respect to the T {R) state). The enlarged inset gives an example of a 
trajectory, in blue, that exhibits multiple passages across the transition region. B. Ensemble 
averages of the RMSD for the T ^ R (top) and R R" (bottom) transitions are obtained over 
50 trajectories. The solid lines arc exponential fits to RMSD/R and RMSD / R" relaxation 
kinetics. C. Time-dependent changes in the angles (measured with respect to the T state) that 
F, M helices make during the T ^ R ^ R" transitions. The inset shows the dispersion of 
individual trajectories for F- helix with the black hue being the average. D. Time-dependent 
changes in the angles (measured with respect to the T state) that K, L hehces make during the 
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T ^ R ^ R" transitions. The inset on the top shows the structural changes in K, L hehces 
during the T ^ R ^ R" transitions. For clarity, residues 357-360 are displayed in space-filling 
representation in white. The bottom inset shows the dispersion of individual trajectories for the 
K- helix. The black line is the average. In C and Y> 9 = cos~^('u(0) • u{t)). 

Figure 4 : T ^ R GroEL dynamics monitored using of two interacting subunits. Side 
views from outside to the center of the GroEL ring and top views are presented for the T 
(left panel) and R (right panel) states. Few residue pairs are annotated and connected with 
dotted lines. The ensemble average kinetics of a number of salt-bridges and contacts between 
few other residues are shown in the middle panel. Distance changes for a single trajectory for 
few residues are given in Fig. S4 in the Supplementary Information. Fits of the relaxation 
kinetics are: {d{t))R58-E209/^ = 14.9 + 9.6(1 - O.lTe-*/^-!^- - 0.83e-*/«25M^), {d{t)) Dsz-Km7 / k = 
8.5+4.9(l-e-*/ioo-0'^^), (ci(t))p33-ivi53/A = 7.3+4.2e-*/6•3/^^ (c^(t))ffi84-i?386/A = 13.2 + 16.5(1- 
0_49e-V20.8M5 - 0.51e-*/85-8'^"), (rf(t))K285-E386/A = 12.6 + 15.8(1 - 0.42e-*/i9.i^s _ 0.516-*/^^-^^^), 
(c?(i))m97-fi386/A = 11.9 + 9.0(1 - 0.296-*/°-^^'^^ - O.Tle-^/^^-^M^), (ci(i)) ^go-Base/ A = 10.4 + 
9.8(0.78e-*/i2.i;..^0.22e-*/6i-«'^^), ((i(i))E257-i?268/A = 9.7+12. l(l-0.35e-*/26-2MS-0.65e-*/66-^'^^). 
Initially, the dynamics of salt- bridge formation between E257 and K321, R322, K245 show non- 
monotonic behavior. Thus, we did not perform a detailed kinetic analysis for these residues. 

Figure 5 : Dynamics of the R —>■ R" transition using two-subunit SOP model simulations. 
The dynamics along one trajectory are shown in Fig. S4 in the Supplementary Information. 
Intra-subunit salt-bridges (or residue pairs) of interest (D83-K327, E409-R501, P33-N153) are 
plotted on the top panel, and inter-subunit salt-bridges (or residue pairs) of interest (E257-K246, 
E257-R268, E257-K321, E257-R322, 1305-A260) arc plotted on the bottom panel. For emphasis, 
K80-D359 salt-bridge dynamics, that provides a driving force to other residue dynamics, is 
specially plotted on the bottom panel. The quantitative kinetic analysis performed for rupture 
of D83-K327 and formation of K80-E359 salt-bridges show {d{t)) Dsz-Ki2'7 / ^ = 10.4 + 26.9(1 - 
g-t/77.9M.)^ (rf(t))^80_^359/A = 14.1 + 26.4e-*/28■0/^^ 

Figure 6 : Transition state ensembles (TSE). A. TSEs are represented in terms of dis- 
tributions P(gi) where = ^ax{RMSD/xy-'t^n(LMSD/x) ■ Histogram in red gives P(gt) for 
T ^ R (red) and the data in green are for the R — > R!' transitions. For T ^ R, 
X ^ R, min{RMSD/X) = 1.5 A and max{RMSD/X) = 8.0 A. For R R", X = R", 
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min(RMSD/X) = 1.5 A and max(RMSD/X) = 14.0 A. To satisfy conservation of the number 
of molecules the distributions are normalized using / dq^ [P{q^\T — > i?) + P{q^\R — > i?")] = 1. 
Twenty overlapped TSE structures for the two transitions arc displayed. In the bottom panel, 
the distributions of txs that satisfy 6^ < 0.2 A, are plotted for the T ^ R and the R R" tran- 
sitions. B. For the T ^ i? TSE we show the salt-bridge distances {d^f^'^^^^, d^l^'^^^^) with 
black dots. The red and the green dots are the equihbrium distances {{d^f^'^^^^) , 
in the T and the R states, respectively. The distance distributions for the TSE are shown in 
blue. 
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SUPPLEMENTARY INFORMATION 



Rationale for using the Self-organized Polymer (SOP) model for allosteric tran- 
sitions: Several studies have shown that the gross features of protein folding l|, |2l], mechanical 
unfolding of proteins and RNA 3|, |J] and the complicated global motions of large systems can 
be captured using native structure-based models |5i] • More recently, such native structure-based 



methods 



10| have also been used to probe transitions between two given end-point 



structures and to ascertain the nature of robust modes that mediate the allosteric transitions 
In order to realistically simulate transitions between distinct allosteric states, it is nec- 
essary to use simple coarse-grained models that can be used to capture, at least qualitatively, 
the inherent dynamics connecting two or more states. Although simple elastic network models 
have given insights into the nature of low frequency dynamics of a number of systems ^ , the 
inherent linearity of the model limits their scope when dealing with potential non-linear mo- 
tions that must be involved in the allosteric transitions 6[ . We have recently introduced a novel 
class of versatile structure-based models that incorporates the fundamental polymeric aspects 
of proteins and RNA. The SOP model, which is easily adopted for use in very large systems, 
has already proven to be successful in obtaining a number of totally unanticipated results for 
forced-unfolding and force-quench refolding of RNA and proteins 

Building on the successful application to single molecule force spectroscopy of biomolecules 
we used a variant of the SOP model to probe the complex allosteric transitions in a prototypical 
biological nanomachine, namely, the E. Coli. chaperonin GroEL. Structures of GroEL (T, i?, 
and R"), that are populated in the reaction cycle, have been determined (see Fig. SI for a side 
view). Given the two allosteric states (say T and R of GroEL) we induce transitions between 
the two states by switching the energy function representing one structure to another. A few 



methods for achieving the switch smoothly have been recently proposed [7|, |8|, |9|, [lO|. We have 
advocated a Langevin dynamics based method in which the switch is accomplished using Eq. 2 
of the main text (see Fig. S2 for a pictorial view). In writing the equations of motions given 
by Eq. 2 in the main text, we assume that initially the ensemble of conformations are in the 
T state so that the conformations obey the Boltzmann distribution i.e, P({rj}) ~ e~^^^^^'^^^'^'> 
where H{{ri\\T) is the SOP energy for the T-state GroEL (see Eq. 2 in the main text). 
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Upon switching the Hamiltonian (see below) over a time interval Nth* the dynamics follows the 
equations of motion given in step (iii) of Eq. 2 in the main text. As long as fluctuation dissipation 



theorem is satisfied, so that the standard potential conditions are met 14|, it can be rigorously 
shown that at long times the biomolecule will reach the R state with the conformations given 
by the equilibrium Boltzmann distribution governed by the Hamiltonian H{{fi}\R). From this 
perspective the procedure we have used is rigorous. 

The assumption that switch occurs over a predetermined duration requires a few comments. 
The allosteric transitions occur as a result of ligand-binding or interactions with other 
biomolecules. As a result of binding, local strain is induced at the interaction site or sites. As 
long as the rate of strain propagation is larger than the rate of conformational change of the 
molecule then the switch over a reasonable time is justified. Extremely rapid switching, which 
is tantamount to very large local loading rates, is unphysical as is adiabatic change in the 
energy function. Given these extreme situations we choose a value of Nt that falls in between 
the extreme conditions. In our procedure Nt can be adjusted to mimic the potential rate of 
strain propagation, which induces the allosteric transitions. The efficacy of the procedure has 
been demonstrated by successful applications to describe the multiple allosteric transitions in 
GroEL. 



Implementation of the Hamiltonian switch: Here we give details of the algorithm for 
executing the second step in the equations of motion (See subsection Additional details in the 
Methods section in the main text). During the transition interval we define r° (T R) using 
linear combination of r°j{T) and r°j{R) where rfj{X) is the distance between the residues i and 
j in the structure X with X = T or R. The switch between r°-{T) and r°-{R) is carried out 
slowly every 100 time steps. In the initial stage of the T ^ R transition (0 to 100 time steps) 
we let r°.(T ^ R) = r°^{T). Subsequently, we set r° (T R) = [l - 0.01A;)r°.(T) + 0.01A;r°.(/?) 
where k = 1,2, 3.... 100 is changed every 100 time steps. Thus, the switch in r° (T R) occurs 
over 10,000 time steps. The loading condition can be varied by changing the number of time 
steps used to achieve the switch in the distances between the native contacts. For convenience 
we used a linear combination of r° (T) and r° (i?) during the switch process. More generally, 
one can use non linear combinations, i.e., r°-{T ^ R) = g{t)r°-{T) + (1 — g(t))r°j{R) where g{t) 
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is an arbitrary function (exponential for example). 

At the end of each interval, the new value of r? (T — > R) is substituted into the SOP 
Hamiltonian to compute the forces needed to solve the equations of motion in step (ii) (see 
Eq. 2 of the main text). In the present application, the procedure for using r°^(T R) lasts 
only < 50 ns. As a result, the dynamics of distant pair is not affected. Only the equilibrium 
distances of native pairs, that are already in contact in the T state but lead to instability in 
the intergration of equations of motion due to rapid switching from T to R, are corrected to 
the equilibrium distances at R state (see Fig. S2). 

Time scales and their relevance: The characteristic time scale of the Brownian dynamics 
in the overdamped limit is th = ^^tl where we used the friction coefficient C, — 50t£"^, 
Eh = 2.0kcal/mol and = (i^)i/2 ^ pg for proteins. The simulations are performed at 
T = 300 K [ksT = 0.6kcal/mol). We chose the integration time step h = O.Itl for (i) and (iii) 
while h* — (0.001 — 0.01)Ti;, for (ii). Thus, 10^ integration time steps with h — O.Itl in our 
Brownian dynamics simulations correspond to 50 /is. Because we have used a minimal model for 
GroEL the time scales quoted in the main text should be viewed as lower limit for the various 
processes. The actual time scales are expected to be much longer. The relatively time scales for 
different aspects of the allosteric transitions are likely to be correct. For example, we find that 
the tilt of K and L helices occur four times more slowly than the F and M hehces (see Fig. 3 in 
the main text). Our prediction of factor of four is, in all likelihood, an accurate estimate. During 
simulations we collected the structures every 0.5 /is to analyze the allosteric transition dynamics. 

Analysis of dynamics: To perform a quantitative analysis on the salt bridge or contact pair 
dynamics we averaged over the time traces of all the trajectories. For the contact dynamics 
of two-subunit GroEL we generated N — 50 trajectories in total and computed dynamic 
changes in specific residue pairs using {d{t)) = J2iLi ^ii^) ■ general {d(t)) is fit using 
{d{t)) = {d{t*)) + A(/e-(*-**)/^i + (1 - /)e"(*"**^/'"'), where t* = 50^s, A is the average decrease 
in the contact distance, and Ti, T2 are the relaxation times for the pathways partitioned into / 
and 1 — /. 
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Figure Captions 

Figure SI: The columns from left to right show the side view of GroEL structure in the T, 
R, and R" states. 

Figure S2: Illustration of the procedure to switch SOP Hamiltonian from T to R state. To 
avoid the computational instability caused by instantaneous switch of equilibrium distance, the 
equilibrium distance from T to R state (r° (T) r°-{R)) is gradually switched using a series of 
transient potentials defined with r° (T R) (see Additional details in Methods section). 

Figure S3: TSEs represented in terms of distribution P(At) where = 1/2 x 
\{RMSD/T){tTs) + {RMSD/R){tTs)\ ior T ^ R transition. for R R" transition is 
similarly defined. 

Figure S4: The dynamical changes in the distances between a number of residues in a 
single trajectory during T — > i? and R — > R" are plotted on A and B, respectively. 
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